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Abstract 

In this work we derive equivalence relations between mimetic finite difference 
schemes on simplicial grids and modiffed Nedelec-Raviart-Thomas finite el¬ 
ement methods for model problems in H(curl) and H{dw). This provides 
a simple and transparent way to analyze such mimetic finite difference dis¬ 
cretizations using the well-known results from finite element theory. The 
finite element framework that we develop is also crucial for the design of 
efficient multigrid methods for mimetic finite difference discretizations, since 
it allows us to use canonical inter-grid transfer operators arising from the 
finite element framework. We provide special Local Fourier Analysis and 
numerical results to demonstrate the efficiency of such multigrid methods. 
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1. Introduction 

We consider mimetic finite difference (MFD) methods for problems in 
H(curl) and iL(div) with essential boundary conditions. Such methods are 
designed in order to have natural discrete analogues of conservation (of mass, 
momentum, etc), symmetry and positivity of the operators. They are also 
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structure preserving discretizations, namely, they form discrete de Rham 
complexes. 

Such discretization techniques were started in the School of A. A. Samarskii 
at the Moscow State University, and they have been further developed and 
analyzed by Shashkov [1] and Vabishchevich [2]. Regarding the MFD meth¬ 
ods, our presentation here follows Vabishchevich [2] and his Vector Analysis 
Grid Operators (VAGO) framework for dual simplicial/polyhedral (Delau- 
nay/Voronoi) grids. 

Many authors have contributed to the research in this held, by applying 
the MFD methods successfully to several applications ranging from diffu¬ 
sion [3, 4, 5], magnetic diffusion and electromagnetics [6] to continuum me¬ 
chanics [7] and gas dynamics [8]. We refer to a recent comprehensive review 
paper by Lipnikov, Manzini, and Shashkov [9] and a recent book by Beirao 
da Veiga, Lipnikov, and Manzini [10] on the subject for details and literature 
review. 

We are interested in the MFD discretizations of two (standard) model 
problems in H(curl) and H{dw). We show that the MFD methods can 
be htted in a more or less standard hnite element (FE) framework which 
leads to convergence results and makes the design of efficient and fast solvers 
for the resulting linear systems quite easy. Our approach is somewhat like 
special discrete Hodge operates and, therefore, is related to the generalized 
hnite diherence approach proposed by Bossavit (see e.g. [11] and references 
therein). We point out that, in the classical hnite diherence setting, conver¬ 
gence results exist, as can be seen in [2], but deriving them is by all means 
not an easy task. Moreover, while we provide details on the constructions in 
2D, the equivalence between the MFD methods and the FE methods carries 
over with trivial modihcations to 3D case as well. We have only chosen 2D 
because it makes the exposition much easier to understand. 

Such connections between the MFD schemes and the mixed FE methods 
for dihusion equations with Raviart-Thomas elements have been already es¬ 
tablished, see [12, 13, 14, 15, 16] and references therein. In fact, designing 
hnite element methods on arbitrary grids is a hot topic and we refer to the 
recent works on agglomerated grids [17, 18, 19] and virtual hnite element 
methods [20, 21, 22]. 

Most of the existing works are on approximation, stability and structure 
preserving properties of the MFD discretizations. Developing fast solvers 
for the resulting linear systems is a topic that needs more attention, since 
the design of fast solvers makes the MFD discretizations more practical and 
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efficient. For FE methods, solvers can be built using the agglomeration 
techniques introduced by Lashuk and Vassilevski [17, 18]. Such techniques 
do not apply to the MFD discretizations (even on simplicial grids!) and, to 
the best of our knowledge, such results are not available in the literature. 
We point out though that on rectangular grids for standard hnite difference 
schemes for i7(div) problems, a distributive relaxation based multigrid was 
proposed in [23]. 

As we have pointed out, our goal is to apply classical multigrid and sub¬ 
space correction techniques [24, 25, 26, 27, 28] for the mimetic discretiza¬ 
tions, by first establishing the relation with Nedelec-Raviart-Thomas ele¬ 
ments. Such approach automatically makes efficient methods such as the 
ones developed by Arnold, Falk and Winther [29] and Hiptmair and Xu 
(HX) [30] preconditioners applicable for the MFD methods. 

Regarding the convergence of W— and R—cycle multigrid with a multi¬ 
plicative Schwarz relaxation proposed in [29], we complement the numerical 
results with practical Local Fourier Analysis (LFA) which provides sharp es¬ 
timates of the multigrid convergence rates. We use a variant of LFA that is 
applicable on simplicial grids (see [31]) and compare the convergence rates 
predicted by LFA with the actual convergence rates of IF—cycle and R—cycle 
multigrid. 

The rest of the paper is organized as follows. In Section 2, we describe 
the MFD schemes on simplicial grids. In Section 3 we derive the “modi- 
hed” Nedelec-Raviart-Thomas FE methods and show their equivalence to 
the VAGO MFD schemes. Section 4 dehnes the multigrid components: 
smoothers, and, with the help of the results from Section 3, the canonical 
inter-grid transfer operators. In this section, we also discuss the setup and 
the design of appropriate LFA for edge-based discretizations and Schwarz 
smoothers. The results obtained from the LFA analysis are shown in Sec¬ 
tion 5, together with the convergence rates of the resulting multigrid algo¬ 
rithm. Finally, conclusions are drawn in Section 6. 

2. Mimetic finite difference discretizations on triangnlar grids 

We consider the following two model problems for u in a two dimensional 
simply connected domain D: 

curl rot u -|- Ku = f, in D, 

— grad div u -|- ku = f, in D, 
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( 2 . 1 ) 

( 2 , 2 ) 



subject to essential boundary conditions (vanishing tangential or normal 
components respectively). We also used u and f to denote solutions and 
right hand sides for both problems without distinguish them in different 
equations and spaces explicitly. The corresponding variational forms (used 
in the derivation of the FE scheme) are: Find u G H(curl) and u G Ff(div), 
respectively, such that 

(rot u, rot v) + fi:(u, v) = (f, v), for all v G H(curl), (2.3) 

(div u, div v) + «:(u, v) = (f, v) , for all v G (div). (2.4) 

In 3D we replace rot with a 3-dimensional curl. In the variational form, 
H(curl) and if (div), are the spaces of square integrable vector valued func¬ 
tions which also have square integrable rot (curl in 3D) or div respectively. 
The functions in the spaces H(curl) and if (div) are also assumed to satisfy 

the essential boundary conditions (u x n) = 0 for (2.3) and (u ■ n) = 0 

for (2.4) where n is the unit normal vector outward to dQ. 

2.1. Mimetic finite differences on a pair of dual meshes 

We consider MFD schemes for (2.1) 
and (2.2) discretized on a pair of pri¬ 
mal (Delaunay) simplicial grid and a 
dual (Voronoi) polyhedral grid. The 
vertices of the Delaunay triangula¬ 
tion are {xf and the vertices of 
its dual Voronoi mesh are the cir- 
cumcenter of the Delaunay triangles. 

We denote the Voronoi vertices by 
and note that each such 
vertex corresponds to a Delaunay 
triangle Dfc, for f = 1,..., Ny. In 
Figure 2.1 we have depicted a pair of 
dual meshes and marked the Delau¬ 
nay grid-points by squares and the 
Voronoi grid-points by circles. As 
is typical in the MFD schemes, we figure 2.1: A Delaunay mesh and its 

assume that all triangles in the tri- Voronoi mesh, 

angulation have only acute angles. 

This assumption guarantees that the Voronoi vertices will always be in the 
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interior of the Delaunay triangles. For 3D analogues of this assumption we 
refer to [2], By duality, to a Delaunay grid point xf, there corresponds a 
Voronoi polygon Vi, 

Vi = {x G D I |x- xf I < |x- xj^l, j = 1,.. .,Nd, j ^ 

and we denote the Voronoi edge Vij = dVi fl dVj. 

We next introduce the spaces of mesh functions associated with the dual 
Delaunay/Voronoi grids. In an MFD scheme, the unknowns are the compo¬ 
nents of u parallel to the edges of the Delaunay triangulation and evaluated 
at the middle of these edges. We orient each of the Delaunay edges by the 
unit vector 

eg = * = 1,. .., Nn, J G >V^(*) = {j \ dV, ^ 0, j = 1,..., Nn}, 

which is directed from the node with the smaller index to the node with 
larger index, see Figure 2.2a. This convention dehnes a function for 

every edge (xg,xg), 

h(bi) = nr-eg, (2.5) 

where ng is the unit normal vector outward to dVi. We then denote by 



(a) Orientation of the De¬ 
launay edges. 



(b) Voronoi polygon 
and notation for the 
divergence operator. 


Figure 2.2 


the set of mesh functions approximating the values of u-eg at the mid-points 
of the edges connecting xg and xg; namely, for alH = 1,... and j such 
that (xg,xg) is a Delaunay edge we set 




u ■ e 


(^g), 
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To complete the MFD approximation of equations (2.1)-(2.2) we now 
introduce the spaces of mesh functions associated with the scalar quantities 
rot u and div u. We modify a little bit of the dehnitions given in [2] to 
serve better our purposes, although essentially we do not change anything 
quantitatively. With the vertices of the Delaunay (resp. Voronoi) grid we 
associate the space of piece-wise constant functions, which are constants on 
the polygons of the Voronoi (resp. Delaunay) grid. We set 

Hd = {m(x) I u(x) = wf, for all X G Vi, i = 1,, JV^}, (2.6) 

ffy = I '*^(x) = , for all X G V/., k = 1,..., Ny}. (2.7) 


In short, the functions in Ho are constants on Voronoi cells and the functions 
in Hy are constants on Delaunay cells. We then dehne the discrete divergence 
operator as following: 


(divfe u)f : = 


1 

meas(V) 


ufj measiVij), 


( 2 . 8 ) 


where meas(V) is the area (volume in 3D) of the Voronoi polygon Vi, and, 
meas(Vj), is the length (area in 3D) of the Voronoi edge (face in 3D) which 
is dual (perpendicular) to the Delaunay edge (xf,x^). The relation (2.8) is 
clearly an analogue of the divergence Theorem on Vi, namely. 


1 

meas(V) 



1 

meas(V)) 



u ■ n. 


In a similar fashion we dehne the discrete operators grad;^ ; Ho —)■ H^i, 
rot/i ; Hy —)■ Ho, and curl/j : Hy —H/j 


(grad;, M)g 

(rot/, u)^ = 
(curl;, m)^ 


riiD _ 

:= (grad/, M)(xg) ■ eg = ri{i,j) \ 

hi 

j) ^g iFj + VU: 0 gg + ufi lu 

meas(Dfc) 

7/'^ - 7/^ 

-V{k,m)^^^, 


(2.9) 

( 2 . 10 ) 

( 2 . 11 ) 


where, in 2D, meas(Dfc) is the area of the triangle with vertices xg, xg and 
xg; /g, /g and //g are the lengths of its edges, (for example, /g = \xf — a;g|); 
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Figure 2.3: Notation used in the definition of rot^ and curl/^ operators 

and = |x^ —xj^|. We refer to Figure 2.3a-2.3b for clarifying this notation. 

Finally, The MFD stencils corresponding to the model problems (2.1)- 
(2.2) are shown in Figure 2.4a-2.4b for a uniformly refined triangular grid. 
Note that the stencils match exactly the ones given in [2]. The modifications 



Figure 2.4: Stencils corresponding to a uniformly refined triangular mesh 
{\H\ = meas(l/i)) 

for the 3D variants of the operators above can also be found in [2] 

3. Equivalence between mimetic finite differences and finite ele¬ 
ment methods 

In this section, for both model problems (2.1) and (2.2), we are going 
to introduce suitable FE methods to derive stencils on arbitrary structured 
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triangular grids that match those obtained by the MFD schemes. The FE 
methods that we consider are based on the variational formulations (2.3) 
and (2.4). For dehnitions of the corresponding Hilbert spaces and results 
on existence and uniqueness of solutions to these model problems, we refer 
to [32]. 

3.1. Finite element discretization for (2.3) 

Next, we recall that the mesh functions are dehned as approximations to 
the tangential components of the solution on the Delaunay mesh. There¬ 
fore, if we would like to construct a FE discretization that matches the 
MFD method from the previous section, it is reasonable to use lowest or¬ 
der H(curl)-conforming Nedelec elements [33, 34] which have the 0-th order 
moments on the edges of the Delaunay mesh of the tangential components 
of u as degrees of freedom. 

In order to approximate the variational problem (2.3) by the lowest-order 
Nedelec’s edge elements, we consider vector valued functions whose restric¬ 
tions on every Delaunay triangle Dk are linear in each component and have 
tangential components that are continuous across the element boundaries. 
Namely, we dehne the Nedelec’s FE space as following 


= W e H(curl) I = 


Oi 

02 


+ h 


y 

—X 


,k = l,2,...,Nv}. (3.1) 


The FE approximation of (2.3) is: Find u/j G such that 

(rotu/„rotVft) + ^ (f.v^,), Vv/, G Vf. 


(3.2) 


As is well known, the degrees of freedom (functionals which uniquely 
determine the elements in V^) are chosen to ensure tangential continuity 
between elements and in the lowest order case, the degrees of freedom are 
the 0-th order moments of the tangential component on each edge, i.e. Uij = 

xD 

Uh ■ e^, with dehned as in Section 2. 

The bases dual to these degrees of freedom have one basis function per 
Delaunay edge (xf,x^), ipij = ^(AjVAj — A^VAj), where A* and \j are the 
barycentric coordinates of the Delaunay grid. In the standard fashion, the 
solution of problem (3.2) is written as u/^ = and the vector of coeffi- 

(ti) 

cients = {uij) is a solution to the linear system of equations . 



Here, the stiffness matrix has elements given by 

Nv 

)(*2j2)(*ljl) ~ 

k=l 

We can now compare the stencil (a row in A^) corresponding to the FE 
discretization and the MFD discretization derived in Section 2. The result 
can be seen in Figure 3. la-3, lb and there is obviously no match. In order 


' Dk 


(rotrot^roj, + n 


(3.3) 




Figure 3.1: Stencils for a general triangulation. 


to see the relation between the MFD and the FE stencils, the key is that 
the MFD degree of freedom uF is basically a scaled Nedelec FE degree of 
freedom and since the degrees of freedom are dual to the basis functions, 
we need to scale appropriately the basis in order to have the same entries in 
the FE stencil. More precisely, if we use the midpoint quadrature rule on 
every edge (which is exact for functions in V^) we have 

W = ■ eg)(a;g) /g (^^(x). (3.4) 

(ij) / (id) 

Hence, an appropriate re-scaling of the basis is (pfj := /g (pij. Using the re¬ 
scaled basis functions, the corresponding stencil is still not in agreement with 
the MFD but the resulting rows are proportional to each other. We easily 

remedy this by scaling also the test functions (pij = Summarizing, 

^km 

the MFD matrix A^^ and the FE matrix A^ satisfy 

(3.5) 


9 














For the right hand side f, in order to be consistent with the FE method, 
we use the following approximation 

^km J O 

Then we have = Dib^ and, hnally, obtain the following relation between 
the MFD solution = (ujj) and the Nedelec FE solution given by 

= ^ {DiA^D 2){UF^) = Dib^, UFF = d^^U^. (3.6) 

Clearly, there is a function in corresponding to the MFD solution which 
we can define as following (in the notation of Section 2) 

(id) 

Using the fact that and cp^j = p>ij we obtain that = 

u/j. As a consequence, from the standard error analysis for the Nedelec FE 
methods for sufficiently regular D we automatically have an error estimate 
and stability for the MED discretization, namely, 

||u - U^^llrot < C/l||f||L 2 (f^), (3.7) 

where ||v||i.ot := a/( rot v, rot v) + k{'v, v) and || ■ ||l 2 (q) is the standard 
norm on D. 

We can also use other approximations of f in the MED schemes. This 
implies to use an approximation f in (3.2) instead of f. It is reasonable to 
assume that ||f — f||L 2 (Q) < Ch, and, therefore, we still have the error esti¬ 
mate (3.7) by the standard perturbation argument and triangular inequality. 

Remark 1. Although the matrix A^^ = DiA^D 2 could be non-symmetric in 
the standard Euclidean inner product, it is symmetric in the inner product 
dehned by {D^^D 2 )- This is a crucial observation which plays an important 
role in the design of efficient solvers for the MED discretizations of prob¬ 
lem (2.1). 

3.2. Finite element methods for (2.2) 

The situation with the model problem (2.2) is a bit more involved. To 
obtain the EE discretization that matches the MFD discretization of (2.2) 
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we borrow some ideas from [15, 16, 19] and construct an FE space on the 
polyhedral grid formed by the Voronoi cells. Let us recall that we have 
Nd Voronoi cells. We consider the space of functions whose divergence is 
constant on these cells. 

= {vh e hr(div) I G iL(div)(V), div v/, = const.}. (3.8) 

While the particular behavior of the functions inside the Voronoi cell does 
not matter for our considerations that follow, just to £x the space, we assume 
that the elements of are in the standard Raviart-Thomas space [35, 34] 
on a sub-triangulation of every Voronoi cell as shown in Figure 3.2 with hxed 
values of the Raviart-Thomas degrees of freedom on dVi for i = 1,... ,N£,. 
With this choice, the hnite element approximation of (2.4) is: Find U/^ G 
such that 

{dwuh,dw\h) + i^{uh,^h) = Vv/jGV)^'^. (3.9) 

There is one degree of freedom associated with each Voronoi edge (x)(,x)]]). 
The corresponding basis functions, dual to these degrees of freedom, are 
dehned using a sub-triangulation of the hexagon as shown in Figure 3.2 
(although this is not necessary, see [16]), and then solving an auxiliary £- 
nite element problem in each hexagon. We denote by RT{H) the standard 

H 



Figure 3.2: Splitting of an hexagon into six triangles to compute the Raviart- 
Thomas basis function on the hexagon. 

Raviart-Thomas space corresponding to a sub-triangulation of a Voronoi cell 
H. If we are on a uniformly rehned grid, this is a hexagon, which explains 
why we denoted a generic cell with H, but in general it is just a convex cell. 
We restrict our considerations to a hexagon H, since the other cases are sim¬ 
ilar. The basis function on H corresponding to the edge (x)(,x})) C dH is 
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determined by solving the following constrained minimization problem: Find 
(pkm G RT{H) {dim RT{H) = 12) such that 


{ ||9^fcm||*-tmin, \\ipkm\\* ^ W^kruh'^, dw ^km = ± meas(H) ’ 

r- 

I ^km ' ^{km),(jl)i ^ dH, 

where || ■ ||* can be the L‘^{H)-norm or any equivalent norm on the space 
RT{H). 

As in the FE discretization of (2.1), which we considered above, to match 
the MFD discretization, we need to re-scale both the test and trial functions. 
We set ipl = lkrn‘fkm, and take as new test functions ipkm = m^km- We 
then obtain that 

A-ft. where (3.10) 

Similarly, properly dehning the approximation of f, we have 

^RTjjRT ^ ^RT 

with = Dib^^. Together with (3.10), we have . There¬ 

fore, the discrete MFD solution corresponds to the following function in 

u™(x) = 

(k,m) 


which also satishes = u/j, namely the MFD solution is also the FE 
solution. Applying then some standard arguments for a sufficient regular 
domain D, we have the following error estimate 

||u - U^^lldiv < C/l||f||L 2 (Q), (3.11) 

where ||v||div := a/( div v, div v) -|- k ( v , v). 

We can also use other approximation f of f and, similarly, obtain the error 
estimate (3.11) by standard perturbation argument and triangular inequality. 
Moreover, we note that a remark analogous to Remark 1 is necessary here 
as well. 
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Figure 4.1: Hierarchy of grids to perform the geometric multigrid method, 
and location of the unknowns 

4. Multigrid solvers for mimetic finite differences 

Our aim is to hnd efficient multigrid (MG) methods for the MFD dis¬ 
cretization of the vector problems. In this section, we focus on problem (2.1). 
A geometric multigrid (GMG) method for the MFD discretizations on the 
Voronoi cells (hexagonal grids) is a topic of our ongoing research and will be 
reported in our future work. We point out, however, the FE framework that 
we established in the previous section provides the necessary conditions for 
applying efficient methods using irregular coarsening strategies and algebraic 
multigrid together with auxiliary space methods (see [17, 30]). Such meth¬ 
ods, however, are not suitable for local Fourier analysis. Since one of our 
goals is to perform the LFA, we focus on the GMG method for problem (2.1) 
here. 

We are interested in applying a GMG method on triangular grids gener¬ 
ated by regular rehnement. In this way, we naturally obtain a hierarchy of 
grids, as shown in Figure 4.1. Next we describe the components for the MG 
algorithm, i.e., smoother and inter-grid transfer operators. 

4.1. Multigrid components 

Smoother. We use a multiplicative Schwarz smoother proposed in [29] as 
the relaxation (smoother) in the standard V- and IF-cycle. This relaxation 
simultaneously updates all the unknowns around a vertex of the Delaunay 
grid as shown in Figure 4.2. Overlapping of the unknowns requires non¬ 
standard tools in the local Fourier analysis for this smoother. 

Inter-grid transfer operators. For the transfer of information between two 
consecutive grids of the hierarchy, canonical inter-grid transfer operators 
based on the FE framework are constructed. We remind that in order to 
obtain the Nedelec canonical prolongation, we need to write the coarse-grid 
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Figure 4.2: Unknowns simultaneously updated in the overlapping block 
smoother, and overlapping of the blocks. 


coarsG-grid edges 

ij 


kl 




fine-grid edges 


(a) Canonical Nedelec prolongation 
matrix. 



(b) Weights for the canonical 
Nedelec prolongation opera¬ 
tor. 


Figure 4.3 


basis functions as a linear combination of the fine-grid basis functions, that 
is, 

coF„(^g)v>£, 

fcZeK 

where K denotes the set of fine-grid edges inside the support of (p^. In 
this way, the coefficients in the linear combination are the entries of the 
prolongation matrix , i.e. {P^){k,i),{i,j) = DOFki{(pfj), see Figure 4.3a. 
Moreover, the weights in P^ are displayed in Figure 4.3b. 

Recall that we needed to re-scale the standard Nedelec basis in order 
to obtain the equivalence between the MFD and FE methods. We mod¬ 
ify accordingly the prolongation operator and, from Nedelec canonical pro¬ 
longation P^, we can write everything in terms of the re-scaled fine-grid 
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1/2 



Figure 4.4: Weights for the canonical prolongation associated with the “mod- 
ihed” Nedelec finite element scheme. 


and coarse-grid basis functions to obtain that the prolongation matrix is 
P = D 2 ^h- The resulting prolongation for a general refined triangle 

with angles a and (3 is given in Figure 4.4. For the restriction, by similar 
re-scaling, we choose 

It is worth noting here that R is the adjoint of P in the inner products induced 
by (T>7^F>2,h) and {Di^hD^\) on the hne and coarse grids, respectively, and 
this agrees with the considerations in Remark 1. 

The choice of these inter-grid transfer operators is, in fact, crucial. We 
emphasize the importance of the relation obtained in Section 3 because with¬ 
out this, it is possible, but, by all means not easy to design efficient GMG 
methods for these MFD discretizations. Regarding the coarse-grid MFD op¬ 
erator, with this choice of the inter-grid transfer operators, direct discretiza¬ 
tion on the coarse-grid results in a “Petrov-Galerkin” coarse-grid operator. 
Namely, note that since = [P^YA^P^, we have 


aFD 


= D,^hA^D2,h = D,^h{PY^A^P^D2, 


H 


= A" D2,,Dll P«D2 ,h 


N 


-1 oN , 


R 


= {D,^u{P^fDY)Al^{D-lP^D2,n) 


= PAl^P. 
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4-2. Local Fourier analysis 

We now briefly use local Fourier analysis techniques to assess the con¬ 
vergence of the resulting GMG method. The LFA (or local mode) analysis 
is introduced by A. Brandt in [24] and is a technique based on the Discrete 
Fourier Transform. To perform this analysis one slightly diverts from the 
boundary value problem in hand and considers periodic solutions on an in- 
hnite regular grid. It is also necessary to have a discrete operator dehned 
with a constant coefficient stencil. The boundary conditions are not taken 
into account, that is, we assume that the boundary effect is negligible. This, 
of course, is not true in general, but in many practical situations the LFA 
gives sharp estimates on the MG convergence rates. 

In the framework of LFA, the current approximation to the solution and 
the corresponding error can be represented by formal linear combinations of 
discrete Fourier modes forming the discrete Fourier space. The LFA then 
identihes invariant subspaces in the Fourier space and studies how multigrid 
components act on these subspaces. A detailed explanation of all varieties of 
local Fourier analysis can be found in [26, 36] and on triangular grids, in [31]. 

We now move on to describe the difficulties in performing LFA for the 
GMG components dehned above. The LFA for our case is indeed nonstan¬ 
dard and we need to deal with several issues described below. 

Simplieial grids. Local Fourier analysis has been traditionally performed for 
hnite difference discretizations on structured rectangular grids. This analysis 
was extended to FE discretizations on general structured triangular [37, 38] 
and tetrahedral [39] grids. The key fact for this extension is to consider an 
expression of the Fourier transform in new coordinate systems in space and 
frequency variables and introduce a non-orthogonal unit basis of chosen 
to ht the geometry of the given simplieial mesh. The basis corresponding to 
the frequencies space is taken as its reciprocal basis and with these settings 
the LFA on simplieial grids is not very different from the LFA on the standard 
rectangular grids. 

Edge-based unknowns. As we saw, the unknowns in MFD discretizations of 
problem (2.1) are located at different types of grid-points, and therefore the 
stencils (the rows of the matrix involve not one, but several different 

stencils. The key is to split the infinite grid into several different subgrids 
in such a way that all nodes belonging to a subgrid have the same stencil, 
and to dehne suitable grid-functions playing the role of the Fourier modes 


16 



for such edge-based discretizations. We refer the reader to [40] for a detailed 
description of such analysis. 

Overlapping Schwarz smoothers. Overlapping block smoothers require a spe¬ 
cial LFA strategy. Classical approaches fail for this class of smoothers, be¬ 
cause an overlapping smoother updates some variables more than once, due 
to the overlapping. The main difhcnlty is that in addition to the initial and 
hnal errors, some intermediate errors appear, and this has to be taken into 
acconnt in the analysis. To our knowledge, there are only few papers dealing 
with LFA for overlapping smoothers, and all of them for discretizations on 
rectangular grids (see [41], [42] and [43]). In [31] this analysis is developed for 
FE discretizations on triangular grids and in [44], a general LFA technique 
on simplicial grids for overlapping smoothers is presented. 

5. Numerical results 

In this section, we demonstrate the efficiency of the GMG method for 
the MFD discretizations based on the mnltiplicative Schwarz smoother and 
the inter-grid transfer operators obtained from the modihed Nedelec FE dis¬ 
cretization. We also show a local Fourier analysis for this kind of discretiza¬ 
tions to conhrm the experimental results. 

We hrst consider problem (2.1) with k = 1 on an equilateral triangular 
domain of unit side-length. In Table 1, we display the smoothing factor p 
and the two-grid p 2 g convergence factors predicted by the LFA, together with 
the asymptotic convergence factor computed by using a hF-cycle on a target 
hne-grid obtained after 10 rehnement levels. Since in practice it is worth to 
know if we can use E-cycles instead of hF-cycles, due to the high computa¬ 
tional cost of the latter, we also show the three-grid convergence factors p 3 g 
predicted by LFA for E-cycles, together with the asymptotic convergence 
factors experimentally obtained. These resnlts are presented for different 
nnmbers of smoothing steps u. From the results in Table 1, we observe accn- 
rate predictions of the asymptotic convergence factors. Moreover, an optimal 
behavior of E-cycle is shown, since the obtained E-cycle convergence rates 
are very similar to the lE-cycle convergence rates. We have seen that in this 
case, the convergence factors are independent of how the smoothing steps are 
distributed, and therefore we do not distinguish different distributions in the 
table. Notice that very good convergence factors, below 0.1, are obtained by 
using a E-cycle with only three smoothing steps. 
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VF-cycle 

H-cycle 

z/ 

p'' 

P2g 

w 

Ph 

P3g 

Ph 

1 

0.462 

0.331 

0.330 

0.337 

0.334 

2 

0.214 

0.124 

0.124 

0.133 

0.132 

3 

0.099 

0.070 

0.069 

0.072 

0.071 

4 

0.046 

0.045 

0.045 

0.052 

0.052 


Table 1: Smoothing (/x), two-grid {p 2 g) and three-grid {p^g) LEA convergence 
factors, together with measured hh-cycle and Id-cycle asymptotic convergence 
rates, p'^ and p^, respectively, for an equilateral triangle and different num¬ 
bers of smoothing steps, v. 


Next, in Figure 5.1, we display the history of the MG convergence for 
different hne grids. We use a ld(2, l)-cycle and the stopping criterion is to 
reduce the initial residual by a factor of 10“^°. As is a well-known property 
of the MG methods, we observe the h-independent convergence behavior. 



Figure 5.1: History of the multigrid convergence of a V{2, l)-cycle for differ¬ 
ent numbers of rehnement levels. 

To study the robustness of the proposed method with respect to parame¬ 
ter K, in Table 2, we show the smoothing and three-grid convergence factors, 
for both W- and H-cycle, predicted by LFA for different values of k and by 
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K, 

p 

VF-cycle (p^) 

Id-cycle (p 3 ^J 

1 

0.099 

0.070 

0.072 

o 

1 

to 

0.099 

0.070 

0.072 

10-^ 

0.099 

0.070 

0.072 

lO-*^ 

0.099 

0.070 

0.072 

10“^ 

0.099 

0.070 

0.072 


Table 2: Smoothing (/i) and three-grid LEA convergence factors by using 
VE-cycle (p^) and Id-cycle (p^) for an equilateral triangular grid, with three 
smoothing steps {u = 3), and different values of parameter n. 



Figure 5.2: Stencil of curl-rot operator for an arbitrary triangulation char¬ 
acterized by angles a and (3 by using mimetic hnite differences. 


using three smoothing steps. From the table, it is clear that the results are 
independent of k. 

To show a more general applicability of the method, arbitrary structured 
triangular grids are considered. These grids can be characterized by two 
angles a and /?, and therefore, after simple computations we can obtain the 
stencil corresponding to curl-rot operator in terms of a and /?, as we can see 
in Figure 5.2. 

In this way, a systematic analysis with the LFA tool can be performed 
for a wide range of triangulations. In Figure 5.3, we display the three-grid 
convergence factors predicted by the LFA for a wide range of triangular 
grids characterized by angles a and (3. For these results, a Id-cycle with 
three smoothing steps has been considered. 

From Figure 5.3, we observe a deterioration of the convergence factor 
when a small angle appears in the triangulation. This behavior is typical 
when point-wise smoothers are considered on a grid with anisotropy. It is 
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Figure 5.3: Three-grid convergence factors for V-cycle predicted by LFA for 
a wide range of triangles characterized by angles a and (3. 


possible to improve these convergence factors by using a relaxation parame¬ 
ter u. For example, if an isosceles triangular grid with angles 80°-80°-20° is 
considered, we obtain a factor of = 0.508, but we can improve this result 
to p^g = 0.252 by considering a relaxation parameter uj = 1.35. These opti¬ 
mal parameters can be obtained for other triangulations by using the LFA. 
Of course, other well-known techniques can be applied to overcome the dif- 
hculties arising from the anisotropy of the grid, such as line-type smoothers, 
but this is not the focus of this work. 

6. Conclusions 

In this work, we showed an equivalence between the MFD schemes on 
simplicial grids and some modihed FE methods. This relation has been 
obtained for two model problems in H(curl) and if(div). This connection 
leads to immediate convergence results for the MFD schemes and also allows 
the construction of efficient multilevel methods using the FE framework. 
Based on the LFA, we theoretically predicted and numerically showed the 
robustness and the efficiency of the GMG method in H(curl) for edge-based 
discretizations on simplicial grids. 
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